seminario_fuentes
Cristina Velazquez y Deyanira Escudero Ramos
r Sys.Date()
LibrerĆas empleadas
##
## Adjuntando el paquete: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union##
## Adjuntando el paquete: 'plotly'## The following object is masked from 'package:ggplot2':
##
## last_plot## The following object is masked from 'package:stats':
##
## filter## The following object is masked from 'package:graphics':
##
## layout## Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUEPregunta 1: Consumo de antibiĆ³ticos en EspaƱa
Determinar el consumo de antibiĆ³ticos por comunidades autĆ³nomas, diferenciando entre recetado y no recetado, para identificar tendencias de uso.
Pregunta 2: RelaciĆ³n entre PIB y resistencia a antibiĆ³ticos en la UE
Analizar si existe una relaciĆ³n significativa entre el Producto Interno Bruto de los paĆses de la UniĆ³n Europea y los niveles de resistencia a antibiĆ³ticos.
## Warning: Missing column names filled in: 'X14' [14]##
## āā Column specification āāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāā
## cols(
## `freq,na_item,ppp_cat,geo\TIME_PERIOD` = col_character(),
## `2012` = col_double(),
## `2013` = col_double(),
## `2014` = col_double(),
## `2015` = col_double(),
## `2016` = col_double(),
## `2017` = col_double(),
## `2018` = col_double(),
## `2019` = col_double(),
## `2020` = col_double(),
## `2021` = col_double(),
## `2022` = col_double(),
## `2023` = col_double(),
## X14 = col_logical()
## )## spc_tbl_ [42 Ć 14] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
## $ freq,na_item,ppp_cat,geo\TIME_PERIOD: chr [1:42] "A,VI_PPS_EU27_2020_HAB,GDP,AL" "A,VI_PPS_EU27_2020_HAB,GDP,AT" "A,VI_PPS_EU27_2020_HAB,GDP,BA" "A,VI_PPS_EU27_2020_HAB,GDP,BE" ...
## $ 2012 : num [1:42] 30 133 30 121 47 170 91 84 124 128 ...
## $ 2013 : num [1:42] 29 133 31 121 46 171 84 86 125 130 ...
## $ 2014 : num [1:42] 30 132 30 121 47 171 81 88 127 129 ...
## $ 2015 : num [1:42] 30 131 31 121 48 171 83 89 124 128 ...
## $ 2016 : num [1:42] 30 130 31 120 49 166 88 89 125 128 ...
## $ 2017 : num [1:42] 30 127 31 118 50 160 90 91 124 130 ...
## $ 2018 : num [1:42] 30 128 32 118 52 158 91 92 124 129 ...
## $ 2019 : num [1:42] 30 126 33 118 53 153 93 93 121 126 ...
## $ 2020 : num [1:42] 31 125 33 118 55 155 91 93 123 133 ...
## $ 2021 : num [1:42] 31 122 34 120 57 157 94 92 120 135 ...
## $ 2022 : num [1:42] 34 124 35 120 62 159 94 90 117 136 ...
## $ 2023 : num [1:42] 35 123 35 118 64 154 95 91 115 127 ...
## $ X14 : logi [1:42] NA NA NA NA NA NA ...
## - attr(*, "spec")=
## .. cols(
## .. `freq,na_item,ppp_cat,geo\TIME_PERIOD` = col_character(),
## .. `2012` = col_double(),
## .. `2013` = col_double(),
## .. `2014` = col_double(),
## .. `2015` = col_double(),
## .. `2016` = col_double(),
## .. `2017` = col_double(),
## .. `2018` = col_double(),
## .. `2019` = col_double(),
## .. `2020` = col_double(),
## .. `2021` = col_double(),
## .. `2022` = col_double(),
## .. `2023` = col_double(),
## .. X14 = col_logical()
## .. )# Cambio de nombre columna
colnames(pibPP)[1] <- "pais"
# Nos quedamos con las Ćŗltimas letras
pibPP$pais <- substr(pibPP$pais, nchar(pibPP$pais) - 1, nchar(pibPP$pais))
View(pibPP)
lista_pais <- list("BE", "BG", "CZ", "DK", "DE", "EE", "IE", "EL", "ES", "FR", "HR", "IT", "CY", "LV", "LT", "LU", "HU", "MT", "NL",
"AT", "PL", "PT", "RO", "SI", "SK", "FI", "SE")
# nos quedamos solo con los paĆses de la UE
pib <- pibPP %>% filter(pais %in% unlist(lista_pais))
View(pib)
# quitar la columna nula
pib <- pib[, colSums(is.na(pib)) < nrow(pib)]
pib## # A tibble: 27 Ć 13
## pais `2012` `2013` `2014` `2015` `2016` `2017` `2018` `2019` `2020` `2021`
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 AT 133 133 132 131 130 127 128 126 125 122
## 2 BE 121 121 121 121 120 118 118 118 118 120
## 3 BG 47 46 47 48 49 50 52 53 55 57
## 4 CY 91 84 81 83 88 90 91 93 91 94
## 5 CZ 84 86 88 89 89 91 92 93 93 92
## 6 DE 124 125 127 124 125 124 124 121 123 120
## 7 DK 128 130 129 128 128 130 129 126 133 135
## 8 EE 74 76 78 76 77 79 82 83 85 86
## 9 EL 71 72 72 70 68 67 66 66 62 63
## 10 ES 91 90 90 91 92 93 91 91 83 84
## # ā¹ 17 more rows
## # ā¹ 2 more variables: `2022` <dbl>, `2023` <dbl>pib_2022 <- pib %>% select(pais, `2022`)
View(pib_2022)
# ggplot del pib en el 2022 (primero lo ponemos en descendente)
pib_2022_desc <- pib_2022 %>% arrange(desc(`2022`))
# sustituir las etiquetas de los paĆses
pib_2022_desc <- pib_2022_desc %>%
mutate(pais = case_when(
pais == "SK" ~ "Eslovaquia",
pais == "SI" ~ "Slovenia",
pais == "EE" ~ "Estonia",
pais == "MT" ~ "Malta",
pais == "LV" ~ "Latvia",
pais == "HR" ~ "Croatia",
pais == "EL" ~ "Greece",
pais == "BG" ~ "Bulgaria",
TRUE ~ pais # Mantiene los nombres que no estƔn en la lista
))Tras el tratamiento de estos datos se puede ver de manera grĆ”fica la comparativa entre los paĆses.
grafico_pib <- ggplot(pib_2022_desc, aes(x = reorder(pais, -`2022`), y = `2022`)) +
geom_bar(stat = "identity", fill = "gold") +
labs(x = "PaĆs", y = "Valor en 2022", title = "PIB por PaĆs en 2022") +
theme_minimal() +
theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas para mejor legibilidad
grafico_pib
Una vez observada la primera variable, generamos los datos del valor medio de la resistencia a antimicrobianos en los paĆses europeos para comprobar si existe una relaciĆ³n o, efectivamente no tiene que ver.
Incidencia_enfermedades <- read_delim("INPUT/DATA/ECDC_encuesta_AMR_incidencia_enfermedades.csv",
delim = ",", escape_double = FALSE, trim_ws = TRUE)## Rows: 98334 Columns: 10
## āā Column specification āāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāāā
## Delimiter: ","
## chr (8): HealthTopic, Population, Distribution, Unit, RegionCode, RegionName...
## dbl (2): Time, CategoryIndex
##
## ā¹ Use `spec()` to retrieve the full column specification for this data.
## ā¹ Specify the column types or set `show_col_types = FALSE` to quiet this message.## # A tibble: 98,334 Ć 10
## HealthTopic Population Distribution Unit Time RegionCode RegionName
## <chr> <chr> <chr> <chr> <dbl> <chr> <chr>
## 1 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 AT Austria
## 2 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 AT Austria
## 3 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 AT Austria
## 4 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 AT Austria
## 5 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 BE Belgium
## 6 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 BE Belgium
## 7 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 BE Belgium
## 8 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 BE Belgium
## 9 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 BG Bulgaria
## 10 Antimicrobial resiā¦ Acinetobaā¦ R - resistaā¦ % 2012 BG Bulgaria
## # ā¹ 98,324 more rows
## # ā¹ 3 more variables: CategoryIndex <dbl>, Category <chr>, Value <chr>## HealthTopic Population Distribution Unit
## Length:98334 Length:98334 Length:98334 Length:98334
## Class :character Class :character Class :character Class :character
## Mode :character Mode :character Mode :character Mode :character
##
##
##
## Time RegionCode RegionName CategoryIndex
## Min. :2000 Length:98334 Length:98334 Min. :1.000
## 1st Qu.:2008 Class :character Class :character 1st Qu.:1.000
## Median :2013 Mode :character Mode :character Median :2.000
## Mean :2013 Mean :2.158
## 3rd Qu.:2018 3rd Qu.:3.000
## Max. :2022 Max. :4.000
## Category Value
## Length:98334 Length:98334
## Class :character Class :character
## Mode :character Mode :character
##
##
## View(Incidencia_enfermedades)
incidencia_2022 <- Incidencia_enfermedades %>%
filter(Time == 2022)
incidencia_2022MF <- incidencia_2022 %>%
filter(Category == 'Male' | Category == 'Female')
View(incidencia_2022MF)
incidencia_2022MF$Value <- as.numeric(incidencia_2022MF$Value)## Warning: NAs introducidos por coerciĆ³nmedia_poblacion_genero <- incidencia_2022MF %>%
arrange(RegionCode, Population) %>% # Ordenar los datos
group_by(RegionCode, Population) %>%
mutate(mean_value = (Value + lead(Value)) / 2) %>%
ungroup()
incidencia <- media_poblacion_genero %>%
mutate(grupo = substr(Population, 1, 3)) %>%
select(-Unit, -Category, -CategoryIndex, -Value, -Distribution)%>%
filter(!is.na(mean_value))
# Ver el dataframe
print(incidencia)## # A tibble: 754 Ć 7
## HealthTopic Population Time RegionCode RegionName mean_value grupo
## <chr> <chr> <dbl> <chr> <chr> <dbl> <chr>
## 1 Antimicrobial resistā¦ Acinetobaā¦ 2022 AT Austria 0 Aci
## 2 Antimicrobial resistā¦ Acinetobaā¦ 2022 AT Austria 3.10 Aci
## 3 Antimicrobial resistā¦ Acinetobaā¦ 2022 AT Austria 0 Aci
## 4 Antimicrobial resistā¦ Acinetobaā¦ 2022 AT Austria 0.909 Aci
## 5 Antimicrobial resistā¦ Enterococā¦ 2022 AT Austria 3.58 Ent
## 6 Antimicrobial resistā¦ Enterococā¦ 2022 AT Austria 8.53 Ent
## 7 Antimicrobial resistā¦ Enterococā¦ 2022 AT Austria 0 Ent
## 8 Antimicrobial resistā¦ Enterococā¦ 2022 AT Austria 90.1 Ent
## 9 Antimicrobial resistā¦ Enterococā¦ 2022 AT Austria 30.8 Ent
## 10 Antimicrobial resistā¦ Enterococā¦ 2022 AT Austria 2.49 Ent
## # ā¹ 744 more rowsView(incidencia)
media_por_bacteria <- incidencia %>%
mutate(grupo = substr(Population, 1, 3)) %>% # Crear una columna con las primeras letras
arrange(grupo, RegionCode) %>% # Ordenar por grupo y RegionCode
group_by(grupo, RegionCode) %>% # Agrupar por grupo y RegionCode
summarise(media_mean_value = mean(mean_value, na.rm = TRUE)) %>% # Calcular la media para cada combinaciĆ³n de grupo y RegionCode
ungroup()## `summarise()` has grouped output by 'grupo'. You can override using the
## `.groups` argument.View(media_por_bacteria)
# Dataframe final con el valor de media final de todas las bacterias
media_por_pais <- media_por_bacteria %>%
group_by(RegionCode) %>%
summarise(media = mean(media_mean_value, na.rm = TRUE))Gracias a estos datos podemos estudiar las bacterias que generan resistencia de manera individual, y ver su incidencia en los paĆses.
# grƔfico que te dice quƩ bacteria tiene mƔs incidencia
ggplot(media_por_bacteria, aes(x = grupo, y = media_mean_value)) +
geom_boxplot() +
labs(title = "Distribution of Mean Incidence by Bacteria Group",
x = "Bacteria Group", y = "Mean Incidence") +
theme_minimal()
# grĆ”fico que te dice quĆ© paĆses tienen esa media para cada bacteria
# Crear el conjunto de datos interactivo con highlight_key
incidencia_keyed <- highlight_key(media_por_bacteria, ~RegionCode)
# Crear el grƔfico ggplot con el texto configurado para el tooltip
scatter_plot <- ggplot(incidencia_keyed, aes(x = grupo, y = media_mean_value, color = RegionCode, text = paste("PaĆs:", RegionCode))) +
geom_point() +
geom_smooth(method = "lm", se = FALSE) +
labs(title = "Population vs Mean Incidence Value by Region",
x = "Population", y = "Mean Incidence Value") +
theme_classic()
# Convertir el grƔfico ggplot en un grƔfico interactivo de plotly con highlight
interactive_scatter_plot <- ggplotly(scatter_plot, tooltip = "text") %>%
highlight(on = "plotly_click", off = "plotly_doubleclick", color = "red", opacityDim = 0.2)## `geom_smooth()` using formula = 'y ~ x'Pregunta 3: Consumo y resistencia en sectores especĆficos en la UE
Se pretende explorar cĆ³mo el consumo de antibiĆ³ticos en diferentes sectores se relaciona con los niveles de resistencia a antibiĆ³ticos en diferentes paĆses de la UE.
Funcion para sacar los datos JSON Tenemos la carpeta āResistencia_Antibioticos_UEā con 27 archivos .jsonld por lo que hemos creado 3 funciones para obtener los datos. Estos datos se corresponden a la ganaderĆa de la UE.
La funcion obtener_nombre obtiene una lista con los nombres de los archivos que contiene la carpeta
La funciĆ³n obtener_archivo carga cada archivo, accede al enlace necesario y descarga la carpeta.zip (con los datos) de cada archivo y las guarda en una carpeta
La funciĆ³n leer_archivo recorre la carpeta que contiene los archivos.zip los descomprime y carga los datos excel. Hace un dataframe con los datos y los almacena en el global enviroment de cada uno de los paises.
library(httr)
library(tidyverse)
library(jsonlite)
library(readxl)
library(rjson)
obtener_nombre<-function(carpeta){
archivos <- as.list(list.files(path =carpeta))
lista_nombres<-list()
for(i in 1:length(archivos)){
posicion1<-regexpr("_", archivos[[i]])
posicion2<-regexpr("\\.", archivos[[i]])
subcadena<-substr(archivos[[i]], posicion1+1, posicion2-1)
lista_nombres[[i]]<-subcadena
}
return(lista_nombres)
}
obtener_archivo<-function(direccion){
lista_paises<-obtener_nombre(direccion)
lista_enlace<-list()
direccion_archivos<-list()
for(i in lista_paises){
cada_pais<-paste0("AMR_",i,".json")
lista_enlace[i]<-cada_pais
}
for(i in lista_enlace){
cada_archivo<-paste0(direccion,"/",i)
direccion_archivos[i]<-cada_archivo
}
for(i in direccion_archivos){
pais<-fromJSON(file= i)
enlace<-pais$links$archive
respuesta_archivo <- GET(enlace)# Hacer la solicitud HTTP para descargar el archivo
nombre_archivo<-basename(enlace)#Extrae el nombre del archivo de la URL
if (status_code(respuesta_archivo) == 200) {# Verificar si la solicitud fue exitosa (cĆ³digo 200, cĆ³digo estĆ”ndar HTTP que significa "OK")
# Guardar el archivo ZIP localmente en formato binario
writeBin(content(respuesta_archivo, "raw"), nombre_archivo)
print("Archivo ZIP descargado correctamente.")
unzip(nombre_archivo, exdir = "carpeta_destino", overwrite = TRUE)
} else {
print(paste("Error al descargar el archivo. CĆ³digo de respuesta:", status_code(respuesta_archivo)))
}
}
}
obtener_archivo("INPUT/DATA/Resistecia_Antibioticos_UE")
leer_archivo <- function(carpeta) {
carpeta_destino <- carpeta
archivos_zip <- list.files(carpeta_destino, pattern = "\\.zip$", full.names = TRUE)
# Iterar sobre los archivos .zip y descomprimirlos
for (archivo in archivos_zip) {
# Descomprimir el archivo .zip
archivos_extraidos <- unzip(archivo, exdir = carpeta_destino, overwrite = TRUE)
print(paste("Descomprimido:", archivo)) # Imprimir cada archivo que se descomprime
# Filtrar el archivo .xlsx entre los extraĆdos
archivo_xlsx <- archivos_extraidos[grepl("\\.xlsx$", archivos_extraidos)]#AquĆ, grepl() busca archivos cuyos nombres terminen con .xlsx
# Verificar si hay algĆŗn archivo .xlsx descomprimido
if (length(archivo_xlsx) > 0) {
# Leer el archivo Excel como dataframe
datos_xlsx <- read_excel(archivo_xlsx[1]) # Leer el primer archivo .xlsx encontrado
# Asignar el dataframe al Global Environment usando el nombre del archivo como variable
nombre_variable <- make.names(basename(archivo_xlsx[1])) # Crear un nombre de variable vƔlido
assign(nombre_variable, datos_xlsx, envir = .GlobalEnv) # Asignar el dataframe al Global Environment
}
}
}
leer_archivo("carpeta_destino")Juntar los dataframes que tenemos en el enviroment 1. Listar los nombres de todos los dataframes que terminan en ā.xlsxā (ajusta si es necesario) 2. Convertir los nombres a una lista de dataframes usando mget() 3. Unir todos los dataframes en uno solo usando bind_rows
eval=FALSE
nombres_dataframes <- ls(pattern = "_AMR_PUB\\.xlsx$")
lista_dataframes <- mget(nombres_dataframes)
df_combinado <- bind_rows(lista_dataframes, .id = "origen")Cargamos el global enviroment En el que tenemos guardado todos los dataframe de cada pais resultante de la funcion leer_archivo, y el dataframe df_combinado del codigo anterior
ModificaciĆ³n y obtencion del dataframe con los datos seleccionamos las columnas que vamos a necesitar y filtramos en la columna de las bacterias solo las patogĆ©nicas.
Cambiamos los nombres de las columnas y lo asignamos al dataframe
paises_UE_df<-df_combinado%>%
select(rep_Country_name,rep_Country_code,zoonosis_name,matrix_name,totUnitsTested,totUnitsPositive,sampUnitType_name,sampType_name,MIC_name,CUTOFFVALUE)%>%
mutate(zoonosis_name = sub(".*", "", zoonosis_name))%>% # Extraer solo la primera palabra
filter(zoonosis_name != "Escherichia coli, non-pathogenic, unspecified")
nuevos_nombres <- c("NombrePais", "Codigo", "zoonosis_name","OrigenMuestra", "TotalMuestras","MuestraPositiva","Tipo_Unidad_Muestra","TipoMuestra","MIC_name","ValorCorte") # Modifica segĆŗn el nĆŗmero de columnas
colnames(paises_UE_df) <- nuevos_nombresMapa Interactivo de positivos en ganaderia por paises de la UE
paises_UE <- c(
"Cyprus", "France", "Lithuania", "Czechia", "Germany",
"Estonia", "Latvia", "Sweden", "Finland", "Luxembourg",
"Belgium", "Spain", "Denmark", "Romania", "Hungary",
"Slovakia", "Poland", "Ireland", "Greece", "Austria",
"Italy", "Netherlands", "Croatia", "Slovenia", "Bulgaria",
"Portugal", "Malta"
)
mapa_mudo <- st_read("INPUT/DATA/mapaMundi") ## Reading layer `ne_10m_admin_0_countries' from data source
## `C:\Users\deyan\fuentes\seminario_fuentes\INPUT\DATA\mapaMundi'
## using driver `ESRI Shapefile'
## Simple feature collection with 258 features and 168 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -180 ymin: -90 xmax: 180 ymax: 83.6341
## Geodetic CRS: WGS 84mapa_mundo_europa <- mapa_mudo %>%
filter(NAME %in% paises_UE)
mapa_mundo_europa <- mapa_mudo %>% filter(NAME %in% paises_UE)
positivos_por_ciudad <- paises_UE_df %>%
group_by(NombrePais) %>%
summarize(total_pruebas = sum(TotalMuestras, na.rm = TRUE),
total_positivos= sum(MuestraPositiva, na.rm = TRUE)) %>%
mutate(ratio_positivo = (total_positivos / total_pruebas) * 100)
# Unir datos de positividad al mapa
mapa_mundo_europa$NAME <- as.character(mapa_mundo_europa$NAME)
positivos_por_ciudad$NombrePais <- as.character(positivos_por_ciudad$NombrePais)
# Realizar el join usando las columnas correctas
mapa_mundo_europa <- mapa_mundo_europa %>%
left_join(positivos_por_ciudad, by = c("NAME" = "NombrePais"))
mapa <- ggplot(mapa_mundo_europa) +
geom_sf(aes(fill = ratio_positivo)) +
scale_fill_viridis_c(option = "plasma", na.value = "gray") +
labs(title = "Tasa de Positividad por PaĆs en Europa",
fill = "Tasa de Positividad (%)") +
coord_sf(xlim = c(-30, 50), ylim = c(35, 72), expand = FALSE) + # Ajustar lĆmites para hacer zoom en Europa
theme_minimal() +
theme(plot.title = element_text(hjust = 0.5))
# Convertir el grƔfico a interactivo con plotly
mapa_interactivo <- ggplotly(mapa)
mapa_interactivo